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The exact form of the Jacobi - Levi-Civita (JLC) equation for geodesic 



' spread is here explicitly worked out at arbitrary dimension for the configura- 



tion space manifold Me = {q G R^V^g) < E} of a standard Hamiltonian 



o 

O ■ system, equipped with the Jacobi (or kinetic energy) metric gj. As the Hamil- 



■ tonian flow corresponds to a geodesic flow on (ME,gj), the JLC equation can 

be used to study the degree of instability of the Hamiltonian flow. It is found 
that the solutions of the JLC equation are closely resembling the solutions 
of the standard tangent dynamics equation which is used to compute Lya- 
punov exponents. Therefore the instability exponents obtained through the 
JLC equation are in perfect quantitative agreement with usual Lyapunov ex- 
ponents. This work completes a previous investigation that was limited only 

to two-degrees of freedom systems. 
PACS numbers(s): 05.45. +b; 02.40.-k; 05.20.-y 
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In recent papers we have investigated the dynamical stability properties of two- 
degrees of freedom Hamiltonians (N = 2) within the framework of a geometric formulation 
of dynamics that makes use of Riemannian geometry. At iV = 2 the phase space structure 
of a system can be investigated in great detail. In fact the use of Poincare surfaces of section 
makes possible to identify the initial conditions that originate regular and chaotic motions 
in the system, so that this qualitative description as well as the measurement of chaos by 
Lyapunov exponents can be thoroughly compared with the outcome of the Riemannian based 
approach. However, the N = 2 case is a very special case, at least from the geometric point of 
view; in fact there is only one curvature function that - at each point - plays the role of scalar 
curvature, Ricci curvature and sectional curvature. Therefore, in absence of any rigorous 
result to extend at arbitrary N the validity of what we found at N = 2, we have explicitly 
studied the large N case and the results are given in the present paper. There is also another 
motivation for the present work. We have recently exploited the Riemannian geometrization 
of newtonian dynamics to analytically compute the largest Lyapunov exponents in large- N 
Hamiltonian systems and, despite some necessary approximation, the analytic results 

are in strikingly good agreement with the numerical results. However, while applying this 
theory to lattice-y? 4 models || we have encountered some difficulties that are now demanding 
adequate improvements. For the sake of simplicity, all the analytic computations were done 
in an enlarged configuration space-time endowed with Eisenhart metric (see below). In this 
framework, the mentioned improvements can hardly be imagined and a richer geometric 
structure, as is the case of (M E ,gj), is needed (we shall better explain why in the sequel). 
Therefore it is of primary importance to check whether the JLC equation on (Mz,gj) fully 
accounts for the degree of chaoticity of the dynamics at arbitrary N. In principle this 
might not be the case: the JLC equation only describes local instability, whereas chaos 
could crucially depend upon some global property of phase space. As a simple example, let 
us think of the Bunimovich stadium (a portion of the plane, bounded by two half-circles 
joined by two parallel lines, where a free particle bounces), where the shape of the boundary, 
being responsible for the mismatch between focusing and defocusing of trajectories, makes 
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the system chaotic. In the case of Hamiltonian flows at N = 2, something similar happens 
when a trajectory approaches the condition V(q) = E: the curvature function becomes 
very large because it contains powers of the quantity [E — V^q)]" 1 and, correspondigly, the 
configuration space trajectories look as if they were reflected by the V(q) = E boundary. At 
large N such a stadium-like effect is no longer present and [E — V(q)] fluctuates around an 
average value with a negligible probability of getting close to zero, therefore "global" effects 
- if any - should work in a subtler way. 

We consider those systems that are decribed by the Lagrangian function (all the indexes 
run from 1 to AT = dimMs) 

L{qA) = \a lk {q)q l q k -V(q) , (1) 
where is the kinetic energy tensor 

a lk (fq k = 2(E -V)=2W. (2) 
Maupertuis' least action principle 

5 f 2Wdt = 5 j {2[E - V{q)]a ik dq i dq k } 1/2 = 5 I ds = (3) 

J 7 J 7 J 7 

variationally defines the natural motions among all the isoenergetic asynchronous paths 7 
joining two fixed endpoints. Hence the arc-length of configuration space is expressed by 
ds 2 = 2[E — V(q)}aikdq l dq k whence g ik — (E — V{q))a ik . In local coordinates the geodesic 
equations on a Riemannian manifold are solutions of the equations 

^ + r,^ = o (4) 

ds 2 3k ds ds y ' 

where s is the proper time and Y l - k are the Christoffel coefficients of the Levi-Civita 
connection associated with g ik . By direct computation, using g ik — (E — V(q))S ik , 
Y l - k = 2^5 tm (djW5 km + d k W5 mj — d m WSj k ) and ds 2 = 2W 2 dt 2 , it can be easily verified 
that the geodesic equations yield 

d£ = "a? ' 1 v (5) 



i.e. Newton's equations derived from the Lagrangian (|l|). These equations of motion can be 
also seen as geodesies of other manifolds J7J besides (Me, Qj)- Among the others, we mention 
a structure, defined by Eisenhart ||, that we have considered with particular emphasis in 
our previous papers |3|-§],|9[] . In this case the ambient space is an enlarged configuration 
space-time M x M 2 , with local coordinates (q°, q 1 , . . . , q N , q 1 ), where (q l , . . . , q ) e M, 
q° E R is the time coordinate, and q N+1 G R is a coordinate closely related to Hamilton 
action; Eisenhart defines a pseudo-Riemannian non-degenerate metric g E on M x M 2 as 

rf S 2 = ^ ^ g ^ = ^ _ 2 y( g ) rfg rfg + rfg ^ dg N+l + dg N+l q dg _ 

Natural motions are now given by the canonical projection 7r of the geodesies of (M x IR 2 , 
on configuration space-time: i : M x I 2 -> M x R. However, among all the geodesies of 
gE the natural motions belong to the subset of those geodesies along which the arclength is 
positive definite 

ds 2 = g^dq^dq" = 2C 2 dt 2 . (7) 

The stability of a geodesic flow is studied by means of the Jacobi - Levi-Civita (JLC) equation 
for geodesic spread. In local coordinates the JLC equation reads as 

where R k ij r are the components of the Riemann-Christoffel curvature tensor. In previous 
papers we have investigated the relationship between geometry and chaos mainly using the 
Eisenhart metric described above. The JLC equation has been used in its exact form with 
Jacobi metric, only in the case of two degrees of freedom systems 0,|T]]: a perfect agreement 
between the description of instability provided by JLC equation and the description of 
instability provided by more conventional methods (Lyapunov exponents, Poincare surfaces 
of section) has been found. Let us now extend our investigation to arbitrary N. To this 
purpose we use a natural chart (in previous works we adopted parallely trasported frames). 
Let us begin by computing the left hand side of Eq.(H). From (VJ k /ds) = dJ k /ds + 
rjy (dq l /ds) Ji we have 
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= -(— + v^A + r*^ (** + v.^jA (9) 

(is 2 ds \ ds v ds J T ds \ ds %3 ds ) 
trivial algebra and the use of Eq.(f|) lead to 

ds* J ~ ds 2 +2L ^ds ds + \ a ' L « + L « L v L » L n) ds ds J ^> 

where di = d/dq\ Then, we use the expression for the components of the Riemann- 
Christoffel tensor to obtain 

and by substituting Eqs.(|ItJ) and QTT]) into Eq.(g) we finally get 

d 2 J k .dtfdP fdT k ri \ dq r dq 



, + 21?,. -f — + hrr -i-^r J1 = ( 12 ) 

ds 2 13 ds ds V dtf J ds ds v ' 

which has general validity independently of the metric of the ambient manifold. Let us now 
derive its explicit form in the case of Jacobi metric. This metric is a conformal deformation 
of the pure kinetic energy metric, i.e. (gj)ij = e~ 2 ^aij. As we are mainly interested in 
studying standard Hamiltonian systems, = 5ij is assumed. For a conformal metric 
{dj)ij = e ~ 2 ^ij one readily obtains the following expression for the Christoffel coefficients: 
T k j = Sjf,i — 5ifj + 5ijf' k , where f yi = dif = df/dq\ Hence Eq.([L^) is transformed into 

and, using the relation ds = e~ 2 ^dt, we can express it in terms of the physical time t instead 
of the proper time s: 

' 2 ( f%, d 4 - f/4) d 4 + ( fM e-V - 2U 4 d 4) P = 0, (14) 



dt 2 V dt dt J dt V dt dt / 

where = df^f . Finally, as the Jacobi metric corresponds to / = ~ln[l/2(i? — V)], it is 

ll ~ 2(E - V) (15) 

2{E -V) 2{E - V) 2 K ' 

e~ 2f = 2(E - V) (17) 



so that the final expression for the JLC equation for (M E ,gj) is 
d 2 J k 1 / dq i _ T ^g fe \ 



„ 2 -^f) * 



1 

+ 



£7- V 



£-V dt dt 



-,k 



J j = . (1* 



Let us now give the explicit form of Eq.(|T2"D in the case of (M x M. 2 ,g E ), the enlarged 
configuration space-time equipped with Eisenhart metric. One easily finds || that only the 
following Christoffel coefficients do not vanish: Tq = {dV/dqi) and = (—dV/dq 1 ), 

hence, using also ds 2 = (dq ) 2 = dt 2 (as we can set 2C 2 = 1), we get 

d.J k d 2 V 

ds 2 dqWq k ^ ^ 

for k,j — l,...,N. The two other components, J° and J N+1 , do not contribute to the norm 
of J and do not enter the evolution equation fll9D , therefore they can be neglected |J. 

It is a very interesting fact that the JLC equation @ yields the usual tangent dynamics 
equation (19) when explicitly worked out for the Eisenhart metric on M x K 2 . On one hand 
we can expect that at least qualitatively Eq . (|T8|) will give similar results to those obtained 
with equation (|ljj), i.e. the usual Lyapunov exponents. On the other hand, the two equations 
(p"8|) and are so different that it is unclear whether a quantitative agreement too has to 
be expected. Geodesies of (M x M?,g E ) project themselves onto geodesies of (M E ,gj): for 
this reason unstable (stable) geodesies of (M x ~R 2 ,g E ) must correspond to unstable (stable) 
geodesies of (Mb, jj). However, no theoretical result guarantees that the average growth- 
rates of the solutions of Eqs. flTSp and flT9p must coincide. We have addressed this point 
by numerically computing the average growth-rates of the solutions of Eqs. fll9f) and (jl8|) - 
let us denote them by Ai and \{ LC respectively - for a given Hamiltonian flow with a large 
number of degrees of freedom; Ai is the conventional largest Lyapunov exponent. 

Numerical computations have been performed for a flow described by the Hamiltonian 

n n . 

np, = £ + £ [ 2 (9i+i - qi)2 + 4 {qi+i ~ qi)4 ] • (20) 

i=l i=l 
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This is the well known Fermi-Pasta- Ulam /5-model [|T0| , a paradigmatic model of non linear 
classical many-body systems extensively studied over the last decades and at the origin of 
remarkable developments in nonlinear dynamics (for instance, the transition between weak 



and strong chaos has been first discovered in this model [O-O] ) 



The numerical integration of the equations of motion (El) derived from the Hamiltonian 



(P0|) has been performed by means of a third order bilateral symplectic algorithm |13fl , and 
the integration of the two stability equations (|19"D and ( |l8l) has been done by means of the 
same bilateral algorithm and of a fourth-order Runge-Kutta scheme respectively. Both Ai(e) 
and \{ LC (e) have been obtained by means of a standard algorithm |14| |, i.e. computing 



A lM ' f> IW*>ll 2 +ll J fa>ll 2 ) (21) 



MM S 1U \,||J(*n-l)|| 2 +||j(*n-l)|| ! 

where t n = nAt, At is some time interval, tjj is the final time such that Ai has attained a 
good "asymptotic" value. 

In Fig.l the values of X{ (e) are compared to the values of Ai(e) and to an analytically 
predicted curve for Ai(e) (see ref. [j|]); e = E/N is the energy density. As the numerical effort 
to integrate Eq.flTSP is heavier than that required to integrate Eq.(|TSp, we computed Ai for 
N = 256 and N = 2000 coupled oscillators, whereas we computed \{ LC for N = 128 and 
N = 256; at A^ = 256 we have only two points that have been computed just as a stability 
check. The excellent agreement between the outcomes of the two stability equations is well 
evident. 

In Fig. 2 the relaxation patterns of Xf LC (t) and of Ai(i) are also displayed. These are very 
similar at high energy density, whereas they show some separation at low energy density: 
the final values are nevertheless always in very good agreement. These results mean that 
equations (|I8"D and ( |T9"D are not - loosely speaking - the "same" equation written in two 
different forms. As a matter of fact, Eq.flT9|) is contained in Eq.(|l8"D so that one could think 
that in some non-trivial way the extra terms cancel out. This is not the case. There are two 
distinct equations to describe the same phenomenon. They are equivalent for what concerns 
the computation of the average instability growth rates of Hamiltonian flows, but they can 



be not equivalent for the further development of the theoretical approach where the average 
curvature properties of the "mechanical" manifolds are linked to the average chaoticity of 
the dynamics through an effective stability equation independent of the dynamics itself [|J. 
In fact Eq. flOf) is valid on (ME,gj), a manifold which has better mathematical properties 
with respect to (M x IR 2 ,^): (Mg, grj) is a proper Riemannian manifold, it is compact, 
all of its geodesies are in one-to-one correspondence with mechanical trajectories, its scalar 
curvature does not identically vanish as is the case of (M x R 2 , g^i), it can be naturally lifted 
to the tangent bundle where the associated geodesic flow on the submanifolds of constant 
energy coincides with the phase space trajectories. 

In conclusion, we have seen that the results found for the N = 2 case j[],0 generalize to 
arbitrary N, hence the phenomenological information given by Lyapunov exponents can be 
retrieved on the manifold (Mg^j) at arbitrary dimension by means of the JLC equation 
for geodesic spread. 
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FIGURES 

FIG. 1. \{ LC (e) computed at N = 128 is represented by full circles and computed at N = 256 
by full triangles. The largest Lyapunov exponent Ai(e) is represented by open circles (iV = 256) 
and open squares (N = 2000). The solid line is the analytic prediction for Ai(e) given in ref. Qj. 

FIG. 2. The relaxation patterns \f LC (t) and \i(t) are compared at different values of the 
energy density. Full symbols denote \± LC (t) and open ones denote Ai(t). From top to bottom 
£ = 392, e = l,£ = 0.075. 
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M.Cerruti — Sola,R.Franzosi,M.Pettini 
'Lyapunov exponents from geodesic...' 

Fig- 1 




M.Cerruti — Sola, R.Franzosi, M.Pettmi 
'Lyapunov exponents from geodesic...' 
Fig. 2 
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